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ABSTRACT 

We present results from the first self-consistent multi-fluid simula- 
tions of chromospheric magnetic reconnection in a weakly ionized reacting 
plasma. We simulate two dimensional magnetic reconnection in a Harris 
current sheet with a numerical model which includes ion-neutral scat- 
tering collisions, ionization, recombination, optically thin radiative loss, 
coUisional heating, and thermal conduction. In the resulting tearing mode 
reconnection the neutral and ion fluids become decoupled upstream from 
the reconnection site, creating an excess of ions in the reconnection region 
and therefore an ionization imbalance. Ion recombination in the reconnec- 
tion region, combined with Alfvenic outflows, quickly removes ions from 
the reconnection site, leading to a fast reconnection rate independent of 
Lundquist number. In addition to allowing fast reconnection, we find 
that these non-equilibria partial ionization effects lead to the onset of the 
nonlinear secondary tearing instability at lower values of the Lundquist 
number than has been found in fully ionized plasmas. These simulations 
provide evidence that magnetic reconnection in the chromosphere could 
be responsible for jet-like transient phenomena such as spicules and chro- 
mospheric jets. 

Subject headings: Magnetic reconnection - Sun: chromosphere - Sun: mag- 
netic fields 
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1. Introduction 



Magnetic reconnection is the general process by which free magnetic energy 
stored in a plasma can be converted into kinetic and thermal energy by breaking 
the frozen-in constraint which exists for a perfectly conducting plasma. Magnetic 
reconnection occurs in a variety of astrophysical plasmas, including the interstellar 



medium (ISM), galactic disks, and the solar atmosphere (Zweibel & Yamada 2009) 



It has been suggested as a cause of many transient phenomena, such as solar flares 



and X-ray jets (Moore et al. 2011), magnetospheric substorms (Borg et al. 2007; Mi 



Ian et al. 2007) and solar 7-ray bursts (Emslie et al. 2005 Tanuma & Shibata 2005) 



Reconnection also occurs in laboratory plasmas (Gray et al. 2010), and plays a key 



role in the self-organization of fusion plasmas (Park et al. 2006; Gangadhara et al. 



2007). 



Magnetic reconnection has long been considered as a mechanism for creating 
transient phenomena in the corona, such as solar flares, coronal mass ejections, and 



more recently X-ray jets (see review by Moore et al. (2011)). In addition, recent 



observations of the solar chromosphere, the cooler, weakly ionized plasma below the 
solar corona, have suggested that reconnection is occurring there also. In particular, 
the quiet chromosphere exhibits localized, transient outflows on a number of different 



length-scales. The largest class of such outflows are called "chromospheric jets" (Shi 



bata et al.|[2007| . These are surges of plasma, observed in Ha and Ca-II, with typical 
lifetimes of 200-1000 s, lengths of 5 Mm, and velocities at their base of 10 km/s. Anal- 
ysis of observations of chromospheric jets with Hinode have shown "blobs" of plasma 
within the jet outflow (Shibata et al. 2007). A smaller class of these outflows are 



called "spicules" (Sterling 2000), and are mainly observed on the solar limb, though 



possible disk counterparts have recently been identified in blue-shifts of Ca I and Ha 



emission (Sekse et al. 2012). Spicules have lifetimes of 10-600 s, lengths of up to 1 



Mm and velocities of 20-150 km/s. 

A unified model of chromospheric outflow generation has recently been suggested 



by Shibata et al. (2007) and Moore et al. (2011). The basic structure of the theoreti 



cal model consists of the so-called "anemone" structure: a bipole field emerging into 
and reconnecting with a unipolar field. By changing the size of the emerging bipole, 
this simple model has been presented as a way of explaining spicules and chromo- 
spheric jets, as well as solar X-ray jets. Other evidence that magnetic reconnection 
is a possible driver of chromospheric transient phenomena includes Alfvenic flows in 
spicules and "blobs" of plasma in chromospheric jets, both of which are by-products 
of magnetic reconnection. 
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A naive estimate of the reconnection rate M (the rate at which magnetic field 
reconnects and ejects plasma) for these chromospheric phenomena can be made by 
dividing the infiow rate of plasma into the reconnection site by the outfiow rate. The 
simple "anemone" structure described above has a fiux infiow reservoir which is of 
the same size as the outfiow region, and when this reservoir of fiux is depleted, the 
reconnection stops. Therefore a typical size L is used for both infiow and outfiow 
regions, and given an outfiow speed Vout^ and a lifetime tuje, the reconnection rate 
can be estimated as M — -, — . Taking the range of lifetimes, lengths and velocities 
for chromospheric jets gives a minimum of M 0.5. Doing the same for spicules gives 
a minimum of M 0.01. 

Traditionally, magnetic reconnection has been considered within the single-fiuid, 
fully ionized, magnetohydrodynamic (MHD) framework, which is applicable to col- 
lisional plasmas. For a highly conducting plasma, the diffusion due to electron-ion 
collisions is negligible and the magnetic fieldlines are frozen into the plasma. Recon- 
nection occurs when the frozen-in constraint is broken on timescales much shorter 
than the classical diffusion time, by allowing fieldlines to reconnect through a narrow 



Parker 


(1957 


) and 


Sweet 


(1958 



reconnection as a local process by considering a current layer of width 6 much smaller 
than its length L, and with non-zero resistivity (77) due to electron- ion collisions. The 
model assumes steady state, i.e., the rate of plasma fiowing into the diffusion region is 
equal to the rate of plasma fiowing out, and takes the length L to be the characteristic 
system length scale. Using Ohm's law for a fully ionized, single fiuid plasma 

E + vxB = r/j, (1) 

and using the plasma momentum equation, a simple scaling law for the reconnection 
rate can be derived: 

M^^^./^^ = ^. (2) 

Here S = ^oVaL/t] is the Lundquist number, and va = B / ^p^iQ is a typical Alfven 
velocity, with B evaluated upstream from the current sheet and p evaluated in the 
current sheet. From this analysis, it can also be shown that the current sheet aspect 
ratio cr = f scales as 

a (X —=. (3) 



Using the one dimensional (ID) semi-implicit model for the quiet Sun of Ver- 
nazza et al. (1981 ) gives a mass density of 7.4 x 10~® kg/m^ for the chromosphere at 1 



Mm above the solar surface. Using a magnetic field strength of 50 G gives an Alfven 
speed of approximately 15 km/s. The temperature at 1 Mm in this ID model is ~ 
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6000 K and the electron density is ~ 10^^ m~^, which gives the electron-ion collision 
frequency to be 10^ s^^, and the Spitzer resistivity to be 0.004 flm. Assuming a typ- 
ical length scale of 0.1 Mm gives a Lundquist number of ~ 10^ and a Sweet-Parker 
reconnection rate of M ^ 10~^. This is slower than the M > 0.01 required if recon- 
nection is to explain the observed lifetimes of both chromospheric jets and spicules. 
The problem is greater in the corona, where S ~ 10^, and the Sweet-Parker model 
predicts reconnection rates much too slow relative to those implied from observations. 
This is the general problem of the Sweet-Parker model: the predicted reconnection 
rate is significantly slower than almost all atmospheric transient phenomena. 



Petschek (1964) modified the Sweet-Parker reconnection formulation by allowing 
plasma to be redirected by standing shock waves which are setup at the ends of the 
diffusion region. This allowed for a shorter diffusion region, with length L', and the 
reconnection rate increased by y^L/L'. The maximum reconnection rate was found to 
For a value of S* = 10^, this gives 0.065: substantially faster than the 



be M 



Sweet-Parker prediction of 0.001, and within the range observed for chromospheric 
spicules, though still too small for chromospheric jets. 

While the Petschek model predicts faster reconnection rates than the Sweet- 
Parker prediction, numerical simulations can only reproduce the Petschek reconnec- 
tion regime if the resistivity is localized near the X-point - so-called anomalous resis- 



tivity due to turbulence in the current layer or ion-cyclotron wave effects (Uzdensky 



2003 


Shay et al. 


2004) 


(Birn et al. 


2001 


. 



More recently, Loureiro et al. (2007) and Huang & Bhattacharjee (2010) con 



sidered the secondary tearing instability of the Sweet-Parker current sheet in a fully 
ionized plasma as an alternative means for obtaining fast reconnection. When the 
current sheet aspect ratio (a = 6/L) reached a critical value of ctc = 1/200, the sheet 
became unstable to the secondary tearing instability. Here, thinner current sheets 
were created between the primary plasmoids, which were themselves then prone 
to breaking into secondary plasmoids. With Sweet-Parker scaling in the laminar 
regime, they found that secondary onset occurred for a critical Lundquist number Sc 
of = 1/cr^ = 4 X 10"^. Above this value, the reconnection rate was independent 
of 5*: M = l/^J~Sc. For a chromospheric Lundquist number of 5* = 10^, the recon- 
nection rate as a result of the plasmoid instability gives M = l/^/S^ = 5 x 10~^, 5 
times the Sweet-Parker prediction. Thus, the secondary instability allows reconnec- 
tion independent of the mechanism which breaks the frozen-in constraint, but still 
gives reconnection rates slower than is needed to explain the parameters associated 
with chromospheric reconnection events. While secondary tearing may be the cause 
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of the plasma "blobs" seen in cliromospheric jets, until now no simulations have been 
performed to see if chromospheric magnetic reconnection could yield these plasmoids. 

This paper studies magnetic reconnection in the partially ionized chromosphere, 
focusing on both fast laminar reconnection and plasmoid formation. The plasma- 
/3, the ratio of plasma pressure to magnetic pressure, can be as high as 10^ and 
as low as 10"^ in this region of the solar atmosphere (Gary 2001). The average 



mass density falls over 4 orders of magnitude in a few Mm (Vernazza et al. 1981) 



from the photosphere to the base of the corona. The ionization fraction, defined 
as the percentage of the plasma which is ionized, ranges from 0.1% to 50%, as the 
neutral density falls off faster than the ionized component density. In the bulk of the 
chromosphere the average collision time between neutral atoms and ions is of the order 



of ms (Vernazza et al. 1981). This is much less than a typical chromospheric time 



scale of 7 s, based on a typical Alfven velocity of 15 km/s and a typical length scale of 
100 km. Hence the ions and neutrals are often treated as a single fluid. However, when 
magnetic diffusion length scales become as small as the neutral-ion collision mean free 
path, e.g., at magnetic reconnection sites, the decoupling of neutrals and ions cannot 
necessarily be neglected, and chromospheric magnetic reconnection should be studied 
in a multi-fluid framework. 

Partial ionization affects magnetic reconnection in several ways. The most obvi- 
ous is the effect on the resistivity. Ion-neutral collisions introduce a Pedersen resis- 
tivity (or equivalently, an ambipolar diffusion) in the single-fluid MHD formulation. 



which acts perpendicular to the magnetic field (Braginskii 1965). Brandenburg & 



Zweibel (1994, 1995) showed that this ambipolar dissipation can create thin current 



structures in weakly ionized plasmas. The ionization level can also affect the recon- 
nection rate through the Alfven speed: for coupled systems, the outflow Alfven speed 
depends on the total plasma density, but for decoupled systems it depends only on 
the ionized component density. 



Previously, Vishniac & Lazarian (1999), Heitsch & Zweibel (2003a), and Lazarian 



et al. (2004) considered magnetic reconnection in the weakly ionized ISM, using a ID 



analytic approach (note that the ID paradigm assumes that outflow is negligible). In 
the reconnection inflow, if the ions that are pulled in by the reconnecting magnetic 
field are decoupled from the neutrals, an excess of ions can build up in the reconnection 
region. Recombination can then act as a sink for the ions, also decoupling them from 
the field. Given a sufficiently large recombination sink, the reconnection rate can 
become independent of the resistivity. This prediction has not been studied in higher 
dimensional magnetic reconnection configurations. Below, we examine this and other 
aspects of magnetic reconnection in a weakly ionized plasma within a two-dimensional 
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(2D) numerical model. 



Smith & Sakai (2008 ) and Sakai & Smith (2009 ) performed two-fluid (ion+neutral) 



simulations of the coalescence of magnetic structures in partially ionized plasmas. 
They found that the reconnection rate for a flxed resistivity decreased as the plasma 
became less ionized, suggesting that jets associated with fast reconnection must occur 
in the upper chromosphere, where the ionization level is higher. However, in their 
investigations, the collision frequency and ionization/recombination rates did not de- 
pend self-consistently on the local plasma parameters (temperature and density), and 
therefore did not vary with space and time during the magnetic reconnection, a phe- 



nomena which is vital to test the predictions of Vishniac & Lazarian (1999), Heitsch 



& Zweibel (2003a), and Lazarian et al. (2004). In addition, their choice of resistivity 



was approximately 10 rim, which is 4 orders of magnitude larger than a typical value 
in the chromosphere, and their Lundquist number was subsequently a relatively low 
value of approximately 10. 

In this paper, the flrst self consistent 2D simulations of chromospheric magnetic 
reconnection in a weakly ionized reacting plasma are performed. The numerical model 
used to simulate this reconnection is presented in §2. The results are presented in 
§3, with particular focus on the scaling of the reconnection rate with resistivity, and 
the onset of the plasmoid instability. In §4 these results are used to evaluate the 
likelihood that magnetic reconnection can explain transient phenomena observed in 
the chromosphere such as spicules and jets. 



2. Numerical Method 

2.1. Multi-Fluid Partially Ionized Plasma Model 

A partially ionized reacting multi-fluid hydrogen plasma model is used to simulate 
reconnection in the solar chromosphere. For a detailed description and derivation of 
the model the reader is directed to Meier (2011) and Meier & Shumlak ( 2012[ ). The 
model is implemented in the implicit, adaptive high order flnite (spectral) element 



code framework, HiFi (Lukin 2008). 



The full model consists of three fluids, ion (z), electron (e), and neutral (n). The 
fluids can undergo recombination, ionization and charge exchange interactions, with 
denoting the reaction rate for interaction r affecting fluid a. 



The recombination reaction rate for ions (the rate of change of ion number density 
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due to recombination), T^'^'^, is defined as 

= -ny^, (4) 

where the recombination frequency 

z/^^" = ne^=2.6 X lO-^Ws"^ (5) 



is obtained from Smirnov (2003) and T* is the electron temperature Tg specified in 
eV. 

The ionization reaction rate for neutrals, r^°", is defined as 

rr = (6) 

where the ionization frequency 

I/-" = n,A—^— ( ^) 6-"^-"/^= m=^s-^ (7) 



-14 



is given by the practical fit from Voronov (1997), using the values A = 2.91 x 10 
K = 0.39, X = 0.232, and the Hydrogen ionization potential 0jon = 13.6eV. Note 
that rf" = -r^°", and T^'^'^ = -F^^^ 

The charge exchange reaction rate, is defined as 

r'^ = ac^{Vcx)ninnVc^ (8) 

where 

V,, = J^vl^ + ^vl^ + vl (9) 
V vr TT 

is the representative speed of the interaction and vf^ = |vj — v„|^ with denoting the 
velocity of species a. The thermal speed of species a is given by vtu = . /^^s^, where 

y 

Tq, is the temperature, rria is the corresponding particle's mass, and ks is Boltzmann's 
constant. Functional forms for the charge exchange cross-section (TcxiYcx) can be 



found in Meier (2011). 



The multi-fluid model can be expressed by taking the neutral continuity equation, 
momentum equation, and energy (or pressure) equation to obtain a set of equations 
for the neutral fluid, and combining the electron and ion versions of these equations 
to obtain a set of equations for the "ionized" fluid. These equations are supplemented 
with Faraday's Law, and the required transport closure equations. Assuming charge 
neutrality in a hydrogen plasma, the electron and ion number densities are set equal 
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to each other {rii = Ue). Also the ionized and neutral atom masses are set equal to the 
proton mass {rrii = rrin = rnp). The resulting system of partial differential equations 
(PDEs) is given below. 

Continuity: 

Due to charge neutrality, only the ion and neutral continuity equations are required. 

fjT) ■ 

^ + v-(n,v,) = rr + rr, (10) 

nil 

^ + v-Kv„) = rr + rr. (11) 

Momentum: 

The electron and ion momentum equations are summed and terms of order (rrie/mpY^^ 
and higher are neglected to give: 



d_ 
di' 

+ r-m,(v„-v,) + R--R-. (12) 



(m^n^v,) + V ■ (m^n^v^v, + P. + Pe) = j x B + Rf + r^^m^v^ - T^'^^m.Vi 



The neutral momentum equation is 
d 

— (min„v„) + V ■ {miUnVnVn + P„) = -Rf + r^m^Vi - r^^iV^ 
ot 

+ r-m,(v,-v„)-R- + R-. (13) 

Here the current density is 

N V X B 

J = eni{-Vi - vj = . (14) 

The momentum transfer RJ^'^ is the transfer of momentum to species a due to identity- 
preserving collisions with species (3: 

^'f = mafiUafapi^ P " V„), (15) 

where map = J^"^^ ■ The collision frequency is given by 



i^ap = npllap\ , (16) 

with Tap = ^^^y^. The cross-section Sj„ = S„j is 1.41 x lO""*^^ m^, and the cross- 
section Sen = ^ne is 1 X 10~^^ m^, assumiug solid sphere elastic collisions (Draine et al. 
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Sje, is assumed to be vrr^ 



1983). The cross section for ion-electron collisions, 
where is the distance of closest approach {e^ / {AneokBTei) with eo the permittivity 
of free space and e the elementary charge). 



The momentum transfer from species /3 to species a due to charge exchange is 



Rex 



Meier 


(2011 


) and 


Meier & 



Shumlak (|2012|), appropriate approximations for these terms are 

4 



pjcx ^ 



-miacx{Vcx) 



vr 



9 
47r 



and 



R^^ ^ miacx{Vcx)ninnVinV^ 



Ti 



IT 



47r 



-1/2 



-1/2 



(17) 



The pressure tensor is Pq, = Pq,I + tTq, where is the scalar pressure and tTq is 
the viscous stress tensor, given by tTq, = — ^^[Vvq + (Vva)"""] where is the isotropic 
dynamic viscosity coefficient for the fluid a. 

Internal Energy: 

Again, combining the electron and ion energy equations together and neglecting terms 
of the order {nie/mpY/'^ and higher gives: 



d 



P. 



at V^*^ 7-1 



V ■ e,Vi + 



7 



Vi • Pi + Ve ■ Pe + hi + h, 



j.E 



+ Vi ■ Rf + Qr - r--miW 



n 2 



ion 
i 



-m. 



{vl - vj) + v„ ■ R- - V, ■ R- + QTn-Q 



cx 
ni ■ 



+ Q 
(19) 



The neutral energy equation is 



dSn 
dt 



+ V ■ {Sn-Vn + V„ ■ P. + h„ 



-v„ ■ Rf + Q 



+ r--m,(v,^ 



v.-r;^ 



(20) 



Here ^/2 + Pa/{l — 1) is the internal energy density of fluid a, and the 

term r*°"0jo„ represents optically thin radiative losses. The ratio of specific heats 
is denoted by 7. Q"^ is the heating of species a due to interaction with species /3, 
which is a combination of frictional heating and a thermal transfer between the two 
populations: = R^^ ■ {yp — Vq.) + maiinai'apiTp — To). The heat fiuxes hg, hj are 
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calculated using the Braginskii closure for a magnetized plasma, and can be written 
as 



Kll^abb + K_L,a(I — bb) 



(21) 



where K\\^a{Ta) and n±,a{na,Ta, |B|) account for the effects of thermal diffusion par- 
allel to and perpendicular to the magnetic field direction (b) respectively, and whose 



functional forms can be found in Braginskii (1965). 



denotes 



The neutral thermal conduction is isotropic h„ = — k„V/cbT, 
thermal energy gain of species a due to a reaction r, with Qf"' = ^ ^ ^ 
Qrec _ Y^'^^ksTi. Q'^a denotes heat flow from species (3 to species a due to charge 



rf^lfcsT™ and 



exchange (Meier 2011 Meier & Shumlak 2012). The temperature of the neutrals is 



given by T„ = Pn/ {nnks)] and the ion and electron temperatures are assumed to be 
equal such that Tj = Pi/ijiikB) = Pe/{nekB] 



T 



Ohm's Law: 

The generalized Ohm's law is given by the electron momentum equation: 
d 



dt 



-ene(E + Ve X B) - V ■ Pe + Rr + Rf - (22) 



In the HiFi implementation all terms on the left hand side of Equation (22) are 
neglected, as in the chromosphere of the Sun, the electron inertial scale c/upe is likely 
much smaller than magnetic diffusion length scales. However, the electron viscous 
stress tensor is preserved in order to represent the effects of microturbulence and 3D 



instabilities (Che et al. 2011), and also to damp the dispersive Whistler and kinetic 



Alfven waves at the shortest resolvable wavelengths. Note that the electron-neutral 



collision term Rg" cannot be neglected. Using the identity Vj 
definition w 



Ve + j/erij, and the 



Vj — Eq. (22) can then be written as 

j X B 



E + (v, X B) = r/j + 



-V ■ P. 



-w, 



en,; 



erii 



where 



(23) 



(24) 



is the electron, or Spitzer, resistivity, which includes electron-ion and electron-neutral 
collisions. The system is closed by the use of Faraday's law ^ = —V x E. 

It is worth noting how the model presented here differs from the partially ionized 



(2006 


) and 


Arber et al. 


(2007) 



fluid approach implicitly assumes that the ions and neutrals are in ionization balance. 
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which is determined by the average density and temperature, and the center of mass 
velocity used in the Ohm's law is an average over ions and neutrals. The Pedersen 
resistivity, present in this partially ionized single-fluid approach, is a consequence of 
this center of mass velocity. The multi-fluid model presented here follows the ions and 
neutrals separately, thus self-consistently including the interactions between ions and 
neutrals which are represented by the Pedersen resistivity in the single-fluid approach. 
A key advantage of the multi-fluid model used in this paper is that ions and electrons 
are allowed to be out of local thermodynamic equilibrium (LTE). This will be shown 
to be vital for the onset of fast reconnection in chromospheric plasmas. 



2.2. Normalization 

The equations are non-dimensionalized by dividing each variable (C) by its 
normalizing value (Cq). The set of equations requires a choice of three normaliz- 
ing values. Normalizing values for the length {Lq = 1 x 10^ m), number density 
(no = 3.3 X 10^^ m~^), and magnetic fleld {Bq = 1 x 10^^ T) are chosen. From these 
values the normahzing values for the velocity (f o = Bq/ ^HQmpTiQ = 1.20 x 10^ m/ s~^), 
time (to = Lq/vq = 0.83 s), temperature (Tq = nipB^/kBHompriQ = 1.75 x 10® K), 
pressure (Pq = B^/fiQ = 0.80 Pa) and resistivity (r^o = HqLqVq = 1.5 x 10^ Qm), can 
be derived. 



2.3. Initial conditions and simplified equations 

The simulation domain extends from -36Lo to 36Lo in the x direction and -QLq 
to QLq in the y direction, with a periodic boundary condition in the x-direction and 
perfectly-conducting boundary conditions in the y-direction. The size of the domain 
has been chosen so that the boundaries and the particular boundary conditions do 
not affect any properties of the reconnection region centered and localized near x = 
0, y = 0. 

The initial neutral fluid number density is 200no (6.6 x 10^^ m~^), and the ion 
fluid number density is no (3.3 x 10^® m~^). This gives a total (ion-|-neutral) mass 
density of 1.11 x 10~^ kg.m"'^, and an initial ionization level {ipi = ni/{ni + of 
0.5 %. The initial electron, ion and neutral temperatures are set to 0.005To (8750K). 
These initial conditions are consistent with lower to middle chromospheric conditions. 



based on ID semi-empirical models of the quiet Sun (Vernazza et al. 1981). 
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In this plasma parameter regime, ion-neutral identity preserving collisions and 
charge exchange (CX) interactions are equally important and have a very similar effect 
of coUisionally coupling neutral and ion fluids with a neutral-ion collision mean free 
path of Lni = VT,n/i^n,i = 140 m. However, the detailed CX physics is substantially 
more complicated and so for simplicity of interpretation CX interactions are neglected 
in this initial study and will be considered in future work. Charge exchange terms 



(terms with the superscript cx) in Equations (12)-(20) are thus dropped 



The ion inertial scale for these plasma parameters is c/upi ~ Im, which is much 
smaller than any scale of interest. Consequently, we neglect the Hall (j x B) and the 



electron pressure tensor (V ■ Pg) terms in Equation (23), the electron viscous stress 



tensor (V ■ vrg) in Equation (12) and Equation (19), and set electron velocity Vg equal 



to ion velocity Vj in Equation (19). Similarly, Malyshkin & Zweibel (2011) showed 



that in the geometry of a reconnection current sheet electron-neutral collisions are 
only important in calculating resistivity, and so the term "^°'^°" w in Equation (23) 



is also dropped. With these simplifications, we solve the following set of governing 
PDEs: 

Continuity: 



dt 
drin 
dt 



+ V ■ iuiVi 



+ V ■ {n„\r, 



-pion I -prec 



-prec _|_ TM01 
n ' n 



(25) 
(26) 



Momentum: 



d_ 
dt 



m^n^v,) + V ■ (m,n,ViV, + P. + Pe) = j x B + Rf + V^m.^Jn - r;^^miv(?7) 



d_ 
dt 



'miUnVn) + V ■ {m.inn-Vn^n + Pn) 



Rm I -prec T^o^m 



(2J 



Internal Energy: 
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+ + V- (^£,v, + ^^ + v,-p, + h, + he^ (29) 

= j . E + V, • Rr + gr - K^'lrn^v^ - QT (30) 
+ rr{\m,vl-cj>,^) + Qr. (31) 



dt 



+ V • (£„v,v, • p„ + h,) = -v„ ■ R^" + 

- rr^m,v^-gr + rr^m,v,2 + gr. (32) 



Ohm's Law: 



E + (v, xB)=77j. (33) 

Note that to investigate the scahng of the reconnection rate with the Lundquist 
number (S), the Spitzer resistivity 77 is made a parameter of the simulations. The 
values used for 77 are [0.5, 1, 2, 4, 8, 20] x 10~^77o, and so 77 hes in the range [0.08,3] Qm. 

The viscosity coefficients for neutrals and ions are set to ii = in = 10~^Co and the 
neutral thermal conduction coefficient is = 4 x 1Q~^kq. The normalizing constants 
are ^0 = fapUoLoVQ and kq — mpnoLoVQ/To. Note that for these plasma parameters 
the isotropic neutral heat conduction is much faster than any of the anisotropic heat 
conduction tensor components for either ions or electrons, with coUisional ion-neutral 
heat exchange g™ ^ V- (hj + hg). Thus, thermal diffusion for all species is dominated 
by neutral heat conduction and is primarily isotropic. 

A Harris current sheet is used for the initial magnetic configuration, and is given 
here in terms of the in-plane magnetic flux A^: 

= --Bo In cosh {y/ A^) (34) 

where B = V x A^e^ and A^ = 0.5Lo is the initial width of the current sheet. 

To provide the outward force to maintain this current sheet, both the ionized 
pressure Pp = Pi + Pe and the neutral pressure P„ are increased in the sheet: 



Pn{y) = ^n + ;^7T77T777^' (36) 



1 1-F 

2 cosh^ (y/Av-) 

where F = O.OIPq is chosen to maintain an approximate ionization balance of 0.5% 
inside the current sheet. These perturbations ensure that the total (ion and neutral) 
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pressure perturbation balances the Lorentz force of the magnetic field in the current 
sheet: 

= -V(P„ + Pp)+j X B. (37) 

At the same time, a relative velocity between ions and neutrals is required, so that 
the frictional force R*" can couple the ionized and neutral fluids and keep the fluids 
individually in approximate force balance: 

= -VPp+j X B + Rf, and (38) 
= -VP„-Rr. (39) 

To give this balance, we choose an initial ion velocity 

(y) = tanh(y/A^) nlvl 

niTinVin cosh^ (zz/A^) Po 

To this initial steady state a small, local perturbation of the flux, Al is added to 
initiate the primary tearing instability and start the magnetic reconnection: 



Al = -ee V*vy e V^W (41) 

where e = O.OIPq-^o- The initial state is shown in Figure [T] Only a subset of the 
domain is shown, and the y axis is stretched by a factor of 40, so that thin structures 
can be clearly seen in the plots. Using the symmetry of the initial conditions, only the 
top right quadrant of the domain is simulated, with the use of appropriate boundary 
conditions. 



3. Results 

3.1. Decoupling of inflow during magnetic reconnection 

As mentioned in the previous section, rj is made a parameter of the simulations, 
and is not dependent on the local plasma parameters, so that a general scaling law of 
reconnection rate with resistivity can be derived. Six simulations are performed which 
have values of the resistivity of [0.5, 1, 2, 4, 8, 20] x 10~^riQ. Firstly, the generic processes 
that are evident in the simulations are highlighted by focusing on one particular 
value of 7] = 0.5 x 10~^t]q. Later, a relationship between reconnection rate and rj is 
determined for the range of rj in these simulations. 

Figure |2] shows the early evolution of the reconnection region resulting from 
the initial perturbation. The Harris current sheet undergoes the tearing instability. 
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magnetic field reconnects at the X-point and fiux is ejected. By t = 537. 5to a Sweet- 
Parker-like reconnection region has been formed. At this stage, ions are pulled in by 
the magnetic field and drag the neutrals with them. 

Figure |3] shows the ion and neutral fiow fields, on a smaller subdomain of the 
simulation, centered on the reconnection region. Panel |3^) shows the vertical veloc- 
ities at time t = 537. 5to, with ions on the top left quadrant and neutrals on the top 
right quadrant. Panel [3|d) shows the magnitude of the horizontal velocity, again with 
ions on the top left quadrant and neutrals on the top right quadrant. The bottom 
two quadrants of panels a) and b) show the fiow vectors for ions, on the left, and 
for neutrals, on the right, both scaled to the same magnitude. The vertical velocities 
in panel [3^) show that the ions and neutrals are decoupled, with the ions fiowing 
faster into the reconnection region than the neutrals. The difference between neutral 
and ion velocity is approximately 90% of the ion velocity. The horizontal velocities 
in panel [sjo) show that the ion and neutral outflows (the strong horizontal velocity 
at x=+/-2.5Lo) are coupled. The difference between the ion and neutral outflow is 
negligible compared to the actual flow. The Alfven speed can be estimated using the 
upstream magnetic fleld of Bup = 0.6i?o, taken at the point in the current sheet where 
the current amplitude reaches half its maximum value, and the total number density 
at the center of the sheet of 280no to give va = 0.035fo. The coupled outflow of ions 
and neutrals, which has a maximum of 0.015fo, is thus approximately half the Alfven 
speed (based on total number density, as the outflow is coupled) at this time. Later 
in time the outflow increases to the Alfven speed. 

This feature of the rj = 5 x lO^^r^o simulation is common to all six simulations. 
The ions and neutrals inflows are decoupled, as ions are pulled in by reconnecting 
magnetic fleld and the neutrals are dragged in via collisions. The timescale of inflow 
is fast enough that the collisions cannot keep the neutrals completely coupled to the 
ions and an excess of ions builds up in the reconnection region, creating an ionization 



imbalance. This is the situation considered for astrophysical plasmas by both Vish 



niac & Lazarian (1999), who treated weakly ionized reconnection with an analytic 



approach, and Heitsch & Zweibel (2003a), who calculated analytic and numerical 



solutions of ID steady state models of weakly ionized reconnection. 



3.2. Reconnection rate scaling with resistivity 

As described in the previous section, the plasma in the reconnection region is 
out of ionization balance as the neutrals are largely left behind by the ions. Figure 
4 shows the four components contributing to ^ in the ion continuity equation. The 
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top left quadrant shows the loss due to recombination, the bottom left quadrant shows 
the loss due to outflow (horizontal gradient in horizontal momentum of ions), the top 
right quadrant shows the gain due to inflow (vertical gradient in vertical momentum), 
and the bottom right shows the gain due to ionization. The color scheme is based 
on a log-scale, and shows that within the reconnection region gains due to ionization 
are negligible relative to the other three terms. Looking at the largest values for the 
remaining three terms, the losses due to recombination and outflow are comparable, 
and add up to equal the gain due to inflow. This shows that the reconnection region is 
close to a steady state, with inflow of ions balanced by comparable contributions from 



recombination and outflow. Recall that the ID models of Vishniac & Lazarian ( 1999 ) 



and Heitsch & Zweibel (2003a) assumed that recombination was fast enough to domi- 
nate over the outflow, and thus the horizontal direction (along the current sheet) was 
ignorable. Instead, Figure |4] shows that for the self-consistently created reconnection 
region in this parameter range the recombination and outflow are comparable and so 
2D effects cannot be neglected. 

The nature of the steady state balance in these reconnection simulations has 
important consequences for the scaling of the reconnection rate. In the standard 
Sweet-Parker reconnection scenario, inflow of ions into the reconnection region is 
assumed to be balanced by outflow of ions. This, along with other assumptions, 
gives the standard scaling of the normalized inflow rate (reconnection rate) of M oc 
^Jl] oc 1/v^- For the weakly ionized plasma in these simulations, the situation is 
very different due to the ionization imbalance. 

The standard steady-state argument of the Sweet-Parker model can be modifled 
to include the reacting multi-fluid equations, in order to derive a scaling relationship 
for the reconnection rate. Figure [5] shows a cartoon of the reconnection region in these 
simulations. The reconnection region, inside which the frozen-in constraint is broken 
by resistivity ?7, has width 6 and length L. There is an external ion density of nj^ex* and 
a current sheet ion density of rii^cs- The inflow of ions into the reconnection region 
is Vin, the outflow of ions is Vout, and the upstream magnetic fleld is B^p. Assuming 
that the system is in steady state, i.e., that ^ = 0, and integrating around the 
reconnection region gives 

ni,extVinL = ni,cs{Svout + SLu'''"' - 5Lz/^°"), (42) 

where u^'^'^ and z/*"" are the recombination and ionization frequencies, defined in Equa- 
tions ([5]) and Defining Uout = Vout/L, and equating Ohm's law evaluated inside 
[E = rjj) and outside {E = VmBup) of the reconnection region, 

v.nBu, = m^'^. (43) 

()/iO 
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to eliminate 5, this steady state equation can be rewritten as 



^'^{iyout + iy''''-iy"^). (44) 

For a plasma in ionization balance, u^^'^ — i/*"", and rii^ext ~ Using these 

relationships and the total momentum equation to derive Vout — va (where va — 
, ), recovers the standard Sweet Parker scahng law 



^ Va \J vaL/iq V 5" ^^^^ 

where B^p is the magnetic field upstream of the current sheet, and ptotai is the total 
(ion + neutral) density in the reconnection region. 

In contrast, in these reacting two-fiuid simulations, the system non-linearly and 
self-consistently forms a current layer where plasma is out of ionization balance, where 
the recombination is comparable to the outflow, and where ionization is negligible 
compared to both recombination and outflow. Hence the inflow rate can be approxi- 
mated by 

''^'''^V<,„, + i.-). (46) 



Note that this equation does not by itself indicate how the reconnection rate will 
depend on S, as it is not clear how iy^^^,ni^cs, ^-nd rii^ext depend on S from the 
equations. We will therefore use our series of simulations performed over a range of 
T] to determine the reconnection rate dependence on S. 

For all six simulations the effective Lundquist number is defined by 

V 

where Lgim = 2.75Lo. Note that Lgim is much smaller than the horizontal extent of 
the domain. We define va,o as the Alfven velocity given a field strength of Bq and 
number density 201no, i.e., based on the initial background magnetic field and number 
density. Note that Ssim varies over simulations due to 7] only. The reconnection rate 
is defined by 

Ms^m = (48) 

Here, jmax is the maximum value of the current density, located at {x,y) — {xj,0), 
within the reconnection region, and B^p is evaluated at {xj, Ssim) where 5 sim is the y 
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location on the line x = Xj at which the current density reaches jmax/'^, i-e., half-width 
at half-max of the current sheet. 

For each simulation, the value of Mgim is taken at a time when the length of 
the current sheet, defined by the distance from the original X-point to the location 
of maximum outflow Vout equals Lgim- This happens at a different time for each 
simulation, but in all cases prior to the onset of any secondary instabilities of the 
current sheet. This instantaneous value of Msim is plotted for each simulation against 
the effective Lundquist number Ssim in Figure [6j The dashed line shows the single- 
fluid Sweet-Parker predicted scaling M oc 1/ v^- Over the range of S which these six 
simulations cover, the reconnection rate is only weakly dependent on the Lundquist 
number, with a slow decrease in M^im with increasing Ssim- The reconnection rate 
is much faster than the Sweet-Parker reconnection rate, as a direct result of the 
decoupling of ions from neutrals and the enhanced recombination in the reconnection 
region. As discussed in section 3.3[ plasmoid formation increases the reconnection rate 



for Ssim > Sc- We also note that each of these simulations have been performed with 
the same viscosity coefficient, and it is known that the reconnection rate dependency 
on resistivity weakens when viscous effects in the reconnection layer become important 
dPark et al.lpSi) ). 



Figure [7] shows the scaling of the reconnection region aspect ratio agim = i^sim/Lsim) 
with Ssim for the six simulations. Also shown is the Sweet-Parker scaling of 1/ v^- 
The simulations do not exhibit the a oc 1/ y/S scaling but show an approximate 
1/5* scaling (the power law fit through the simulation points has an exponent of 
— 1.1 ± 0.17). This 1/5 scaling is consistent with the result of Figure [6| that Mgim is 
approximately independent of Ssim- This is shown as follows: Assume that Vout does 
not substantially vary over different simulations, and note that Ssim only changes over 
simulations due to r]. Then using the definition of Mgim and jmax ~ Bup/Ssim gives 

■^sim ^ jmax / -Blip ~ ^/ ^sim OC 1/ Ssim^^sim ~ COUSt. 



3.3. The secondary (plasmoid) instability 

As discussed above, at a critical aspect ratio of ac = 1/200, a resistive current 
sheet can become unstable to a secondary tearing instability known as the plas- 
moid instability. Of the six two-fiuid simulations in this paper, three cases show 
evidence of secondary tearing. Figure [8] shows the onset of the plasmoid instability 
for the simulation with 7] = 0.5 x lO'^rjQ, with the two contours of magnetic fiux 
at A, = -0.06815oLo (white line) and A, = -0.0687SoLo (black line) following 
the evolution of two particular field-lines in time. The initial laminar current sheet 
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breaks up into a number of plasmoids, with thinner current sheets between them. The 
current density and recombination rate within this fragmented reconnection region 
are shown in Fig. |8] on a pseudo-color scale. During this phase of the reconnection 
the formation of a plasmoid chain redistributes the ionized plasma within the cur- 
rent sheet, which is otherwise approximately uniformly distributed throughout the 
reconnection region, into regions of higher and lower electron number density. Since 
r^^'^ oc Ug, the recombination rate is increased within the plasmoids with respect to 
the sub-layers between them. This leads to the plasmoid magnetic flux collapsing 
on itself (i.e, disappearing as it would in vacuum) at a rate comparable to or faster 
than the plasmoids are exhausted out of the reconnection region. This collapse can 
affect the expected distribution of plasmoid sizes and the ionization fraction within 
the plasmoids relative to the background medium in the reconnection exhaust. 



Huang & Bhattacharjee (2010) found that in a fully ionized plasma, the onset of 
the plasmoid instability occurred at a critical value of Lundquist number of 4 x 10^ 
which, with the Sweet-Parker scaling of current sheet width with Lundquist number, 
corresponds to a critical current sheet aspect ratio of 1/200, as shown by the dashed 
line in Figure [7j The three simulations that undergo secondary tearing are shown 
as diamonds in Figure [7j The intersection of a power law fit to the asim{Ssim) data 
intersects the agim = 1/200 line at a value of Sgim = 10^. The simulation with 
Ssim = 10*^ exhibits the plasmoid instability, as do the two simulations with higher 
Ssim (shown as diamonds in the Figure). These two facts support the postulation 



of Loureiro et al. (2007) that the criterion for the plasmoid instability is that the 



aspect ratio decreases below a critical value. For our simulations, as for [Huang 



Bhattacharjee (2010), this critical value is 1/200. 



To demonstrate the change in reconnection rate due to the plasmoid instability, 
Mgim is evaluated again later in each of the three plasmoid-unstable simulations. It 
is difficult to measure Msim during the plasmoid instability, as the reconnection rate 
varies with time, depending on the plasmoid evolution. To give an idea of how the 
plasmoid instability is affecting the reconnection rate we plot in red in Figure [6] the 
range of Msim observed after plasmoid formation in each simulation until the current 
sheet width becomes spatially unresolved. It is apparent, particularly in higher Ssim 
simulations, that the reconnection rate increases as the plasmoids develop. Thus, the 
reconnection rate in this regime is determined by both the fast recombination of ions, 
and the effect of the secondary tearing instability. 
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4. Discussion 



In this paper magnetic reconnection in tlie solar chromospliere is simulated using 
a partially ionized reacting multi-fluid plasma model. The number densities and 
ionization levels are consistent with lower to middle chromospheric conditions. The 
dependence of reconnection rate on Lundquist number is investigated by setting the 
resistivity to be a parameter of the simulations. A simple 2D Harris current sheet 
conflguration with a local perturbation to the in-plane flux is used. The system self- 
consistently and non-linearly creates a reconnection region which is out of ionization 
balance, due to the decoupling of ion and neutral inflows. The model is able to capture 
this physical effect as the two fluids are followed separately, and the ionization and 
recombination rates are self-consistently calculated based on local plasma parameters. 
In the reconnection region, recombination of excess ions is comparable to the outflow 
of ions, which leads to fast reconnection independent of the Lundquist number, and 
when normalized properly, the reconnection rate is approximately 0.1. It is worth 
noting that guide (out of plane) magnetic fleld, not included in these simulations, 
will have an effect on the scaling, as the flux associated with the guide fleld can 



inhibit fast reconnection (Heitsch & Zweibel 2003b) 



The onset of fast reconnection in weakly ionized astrophysical plasmas was pre- 



dicted by Vishniac & Lazarian (1999) and Heitsch & Zweibel (2003a), assuming that 



recombination dominates outflow in the current sheet so the system could be treated 



with one dimensional analysis. In this regime Heitsch &; Zweibel (2003a) found that 
fast reconnection occurred when 



/3q trecombtvi 



10 



/2 



(49) 



where was either 10~^ for the numerical solution, or 10~^ for the theoretical solu- 
tion, (5q = g2/^° ; tree is the tlmescale of recombination, = L'^/rj, and Iad = L'^/Vp- 
In our multi-fluid simulations, the smallest resistivity is ~ 0.08 Qm, rjp ^ 0.3 Qm, 
and tree = l/i^^^^ ~ 30 s, and so Z ^ 4 x 10~^, which lies in the 'fast' reconnection 



regime of Heitsch & Zweibel (2003a). While we have shown that while the prediction 



of fast reconnection by Vishniac & Lazarian (1999) and Heitsch & Zweibel (2003a) 



holds in the chromosphere, in 2D reconnection the recombination is not the domi- 
nant mechanism for the removal of ions from the reconnection region, as outflows are 
equally important. 

The critical aspect ratio at which secondary tearing sets in, ae = 1/200, is reached 
at a lower Lundquist number than has previously been seen in single-fluid simulations 
with /3 ~ 1 (Huang & Bhattacharjee 2010 Samtaney et al. 2009). This is because 
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in this fast reconnection regime a oc l/S*, compared to Sweet-Parker reconnection 
where a oc \l\fS. Note that Ni et aL (2012) found that for /3 = 50, the onset occurs 



at a Lundquist number of 2000-3000, and an aspect ratio of 1/60. The plasmoids 
which form due to secondary tearing in these simulations are losing ions due to re- 
combination, and so their evolution is potentially very different from those seen in 
fully ionized simulations. The further evolution of the plasmoid beyond the initial 
formation and collapse seen here is left to a follow-up investigation. 

It has been conjectured that current sheets are ubiquitous in the chromosphere 



(Goodman & Judge 2012). It has also been established that spicules and chromo- 



spheric jets are ubiquitous solar phenomena (Sterling 2000). Magnetic reconnection 



in chromospheric current sheets is therefore a promising mechanism to link these 
two ubiquitous phenomena, and would explain the formation of spicules and chro- 
mospheric jets. These simulations show that the chromosphere can exhibit fast re- 
connection with rates that are comparable to the estimated reconnection rates of 
observed outflows, as well as the creation of plasmoids, or "blobs" of plasma, due 
to the secondary tearing mode. The reconnection outflow velocities we find in these 
simulations are 5 km/s, which is close to the speeds observed for chromospheric jets, 
but is small compared to the speeds of 20-150 km/s observed for spicules. Using 
the scahng argument that Vout oc B^p/ ^/rii + Un, we expect that by increasing the 
magnetic field from lOG to 50G, and by moving higher up in the chromosphere where 
the neutral density is an order of magnitude lower (and the ionization level increases 
from 0.5% to 10%) spicule-like outflow velocities of up to 80 km/s can realistically be 
achieved. 

The current sheets formed in these simulations of chromospheric reconnection are 
out of ionization balance to such a degree that short recombination times of 30 s are 
created. These times are much smaller than the estimated 3 minute recombination 



time of an acoustic shock heated chromosphere (Garlsson & Stein 2002), highlight 



ing magnetic reconnection as the prime transient phenomena in the chromosphere 
driving non-LTE recombination. The major advantage of the multi-fluid approach 
over single-fluid models is that it can self consistently produce non-LTE high ion 
density structures. This paper has shown that such ion density structures form in 
magnetic reconnection, and that the resulting fast recombination affects the reconnec- 
tion physics. Hence multi-fluid simulations such as these are vital to understanding 
transient phenomena in the chromosphere. 



Finally, it should be noted that Shay et al. ( 2004 ) and others have argued that the 



addition of the ion inertial and/or finite Larmor radius effects to single-fluid resistive 
MHD is essential to obtain the resistivity-independent "universal" fast reconnection 
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rate of M ki 0.1. The simulation results presented here yield this same resistivity- 
independent reconnection rate without any of the ion-electron decoupling effects; 
relying instead on the combination of enhanced recombination rate with generation 
of the secondary plasmoids in the reconnection region. 
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Fig. 1. — Initial conditions for a subset of the whole domain. Note that the y- 
coordinate is expanded by a factor of 40. The top left quadrant shows ionized pressure 
((Pj + Pe)/Po), the top right quadrant shows neutral pressure (Pn/Po), the bottom 
left quadrant shows current density {j / {Bq/ fioLo)), and the bottom right quadrant 
shows vertical ion velocity (vi^y/vo). The solid lines are 10 contour levels of the flux 
Az evenly distributed in the interval [—0.04, — 0.015]So-^o- 
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Fig. 2. — Formation of Sweet-Parker current sheet for the simulation where i] = 
0.5 X 10~^?7o. Current density {j / Bq/ {hqLq)) is shown in the color contours, and 10 
contour lines show A^, evenly distributed in the interval [—0.04, — 0.015]i?o-Z^o- 
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Time=537.5tQ 



Fig. 3. — Plots of ion and neutral flow, showing inflow decoupling and outflow cou- 
pling during reconnection. Panel a) shows color contours of the vertical velocity 
{vi^y/vo) for the ions (top left quadrant) and for the neutrals {Vn^ylvQ, top right quad- 
rant). The contour lines show 10 values of A^,, regularly distributed in the interval 
[—0.04, —Q.Q1\BqLq. The arrows on the bottom left quadrant represent the plasma 
flow, and those on the bottom right quadrant the neutral flow. Panel b) shows the 
same as panel a) but for the magnitude of horizontal velocity for ions (|fi,a;/'yo|) and 
neutrals {\vn,x/vo\). 



- 28 - 



4.0e-02 



1.0e-02 




4.0e-05 



Fig. 4. — The steady state reconnection region showing contributing sources and 
sinks of ions in the current sheet (in units of Lg'^tg ^): Top left quadrant shows rate 
of loss of ions due to recombination. Bottom left quadrant shows rate of loss of ions 



due to outflow 



driiVi, 



Top right shows rate of gain of ions due to inflow 



driiVi 



Bottom right quadrant shows rate of gain of ions due to ionization. The solid lines 
are 10 contour lines of Az, evenly distributed in the interval[— 0.0378, — 0.037]-Bo-^^o- 
This shows that loss of ions due to recombination and outflow are comparable, and 
combine to balance the inflow of ions, with ionization playing an insignificant role. 
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Fig. 5. — Cartoon of the Sweet-Parker reconnection scenario. The current sheet 
length is L and width is 5. The magnetic field just upstream from the current sheet 
is Bup- The ion infiow and outfiow speeds are f and Vouti respectively. The external 
and current sheet ion number densities are rii^ext and rii^cs^ respectively. 
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Fig. 6. — Normalized magnetic reconnection rate M^j^ for 6 simulations with different 
Lundquist number (Ssim)- The squares show the reconnection rate taken at a time in 
each simulation when the length of the current sheet has reached Lgim = 2.75Lo, and 
in all simulations is before the onset of the plasmoid instability. The red lines show the 
range in reconnection rate taken again at later times in the three plasmoid-unstable 
simulations, after the plasmoids are formed. The dashed line is the Sweet-Parker 
scaling law M oc The dot-dashed line shows the separation between plasmoid- 

stable and plasmoid-unstable regimes for these multi-fluid simulations, and is given 
by the value of Ssim where the power law fit of the simulation data <7sim{Ssim) nieets 
the asim = 1/200 line in Figure [Tj 
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Fig. 7. — Scaling of simulation current sheet aspect ratio agim with Lundquist number 

{Ssim)- The squares are simulations where no secondary (plasmoid) instability is seen, 
and the diamonds are simulations where plasmoids are observed. The dotted line is 
the theoretical aspect ratio at which the plasmoid instability sets in. The dashed line 
is the Sweet-Parker scaling law (oc \Jl/ S). The solid line shows a line of best fit of 
the data to a power law. The exponent in the power law is —1.1 ± 0.17. This solid 
line intersects the Usim — 1/200 line at approximately Ssim = 10^. 
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Fig. 8. — Plasmoid formation and evolution: Recombination rate T^^'^Li^tQ (left) and 
current density j^QL^/ Bq (right), and two contour levels of of — 0.0681i?o-^^o (white) 
and — 0.0687i?o-Z^o (black), at three different times in the simulation where t] = 0.5 x 
lO-^r^o. 



